This is a follow-up to my previous post on where Madison is most dense. But what do these dense areas look like? Inspired by the “Every Lot” bot, I’ll use the static Google Street View API via the googleway package to (hopefully) create some visual sense of density.

## Code from previous post to get 10 most dense block groups in Dane County
library(sf)
## Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
library(tidycensus)
library(tidyverse)
## -- Attaching packages ---------------------------------- tidyverse 1.2.1 --
## v ggplot2 3.2.1     v purrr   0.3.3
## v tibble  2.1.3     v dplyr   0.8.3
## v tidyr   1.0.0     v stringr 1.4.0
## v readr   1.3.1     v forcats 0.4.0
## -- Conflicts ------------------------------------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(ggplot2)
library(leaflet)
## Warning: package 'leaflet' was built under R version 3.6.2
library(googleway)
## Warning: package 'googleway' was built under R version 3.6.2
key <- Sys.getenv("STREET_VIEW_API_KEY") #key requires a credit card on file

pop <- get_acs(geography = "block group", 
                     variables = "B01003_001", 
                     state = "WI",
                     county = "Dane",
                     geometry = TRUE,
                     keep_geo_vars = TRUE)
## Getting data from the 2013-2017 5-year ACS
## Downloading feature geometry from the Census website.  To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.
pop2 <- pop %>% 
  mutate(pop_density = estimate / (ALAND /27878400)) %>%  #create variable for population density
  drop_na %>%
  filter(pop_density > 100000) %>% 
  arrange(desc(pop_density)) %>% 
  head(10)

A good starting point for exploration would be the centroid of the block group. Block groups can have odd shapes, though, which can lead to the true centroid ending up outside the polygon. I’ll use st_point_on_surface instead.

pop2 <- pop2 %>% 
  mutate(center = st_point_on_surface(geometry))
## Warning in st_point_on_surface.sfc(geometry): st_point_on_surface may not
## give correct results for longitude/latitude data

A quick visual check:

library(ggmap)
## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
## Please cite ggmap if you use it! See citation("ggmap") for details.
ggplot(data = pop2) +
  geom_sf(aes(fill = pop_density), 
          alpha = 0.5,
          color = NA) +
  geom_sf(aes(geometry = center))

Looks good! Turns out that ten densest block groups form one contiguous area. Now let’s get Street View images. Testing with the first row:

google_streetview(location = st_coordinates(pop2$center[1]),
    size = c(400,400), output = "plot",
    key = key)

“No imagery here” – that’s because lat/lon are in the reverse order! I’ll fix this, and in order to save API calls I’ll also turn on the response_check option in google_streetview

google_streetview(location =  c(st_coordinates(pop2$center[1])[2],
                                st_coordinates(pop2$center[1])[1]),
    size = c(400,400), 
    output = "plot",
    response_check = TRUE,
    key = key)
## Warning in google_streetview(location = c(st_coordinates(pop2$center[1])
## [2], : No imagary was found at location 43.071489,-89.3983367307066

Still no luck. What apparently is going on is that there is no Street View imagery at the exact coordinates – which isn’t too surprising that the coordinates are just some point inside a block group and not necessarily on or near a street. The Every Lot project doesn’t have this issue, as per definition a lot has a street address. A couple ideas to resolve the problem:

The first option is easier to implement and I’ll try it out first. This XKCD comic helps figure out how much rounding is needed:

![XKCD comic on coordinate precision(https://www.explainxkcd.com/wiki/images/8/88/coordinate_precision.png)

“A specific suburban cul-de-sac”-level accuracy with 3 digits sounds good.

google_streetview(location =  round(c(st_coordinates(pop2$center[1])[2],
                                st_coordinates(pop2$center[1])[1]), 3),
    size = c(1000,400), 
    fov = 120,
    output = "plot",
    response_check = TRUE,
    key = key)

Great, it works! I also adjusted the fov (field of view) to 120 degrees in order to get as much context as possible in the image.

Time to pack this into a function and get images for the top-10 densest block groups.

get_sv <- function(x) {
  print(paste0("Block Group ", pop2$BLKGRPCE[x], " in Census Tract ", pop2$TRACTCE[x], " has a population density of ", round(pop2$pop_density[x], 0), " people per square mile and it looks like so:"))
  google_streetview(location =  round(c(st_coordinates(pop2$center[x])[2],
                                st_coordinates(pop2$center[x])[1]), 3),
    size = c(1000,400), 
    fov = 120,
    output = "plot",
    response_check = TRUE,
    key = key)
}

10th-densest block group: Block Group 2, Census Tract 16.06

## [1] "Block Group 2 in Census Tract 001606 has a population density of 326882 people per square mile and it looks like so:"

9th-densest: Block Group 1, Census Tract 16.05

## [1] "Block Group 1 in Census Tract 001605 has a population density of 369497 people per square mile and it looks like so:"

8th-densest: Block Group 2, Census Tract 16.04

## [1] "Block Group 2 in Census Tract 001604 has a population density of 430858 people per square mile and it looks like so:"

7-th densest: Block Group 1, Census Tract 16.04

## [1] "Block Group 1 in Census Tract 001604 has a population density of 495672 people per square mile and it looks like so:"

6th-densest: Block Group 1, Census Tract 16.03

## [1] "Block Group 1 in Census Tract 001603 has a population density of 627955 people per square mile and it looks like so:"

5th-densest: Block Group 3, Census Tract 16.04

## [1] "Block Group 3 in Census Tract 001604 has a population density of 662878 people per square mile and it looks like so:"

4th-densest: Block Group 1, Census Tract 16.06

## [1] "Block Group 1 in Census Tract 001606 has a population density of 688710 people per square mile and it looks like so:"

3rd-densest: Block Group 2, Census Tract 16.03

## [1] "Block Group 2 in Census Tract 001603 has a population density of 749637 people per square mile and it looks like so:"

2nd-densest: Block Group 4, Census Tract 16.04

## [1] "Block Group 4 in Census Tract 001604 has a population density of 907845 people per square mile and it looks like so:"

The densest of them all: Block Group 4, Census Tract 16.06

## [1] "Block Group 4 in Census Tract 001606 has a population density of 1307631 people per square mile and it looks like so:"

Conclusion

Static Street View images of a single location in a block group don’t tell the whole story of density. But they do provide a nice snapshot of the built form in the dense parts of Madison: Pretty much all of the buildings you can see cater primarily to UW-Madison students, but they take on different forms. You have freestanding historic multiplexes with probably 3-10 units (BG 1, tract 16.05) that may once have been large single-family homes; generic 3-5 story apartment buildings from the second half of the 20th century; and then mid- and high-rise buildings, either built as residence halls by the university (e.g. Witte Hall with over 1200 units in the densest block group) or as private development (e.g. Lucky Apartments in the third-densest block group). Block group 1 in tract 16.06 does a great job of showing all of these building types in one shot.